!-------------------------------------------------------------------------------------
subroutine set_ic( dfu, uxMat, uyMat, f0vec, finit, it)
!-------------------------------------------------------------------------------------

  use dimensions
  use mpi_wrap

  use fft_3d

  use parameters, only : kz, vPhase, ux0, uy0, &                      ! real,     phys
                         nz, iWaveFrame, iModx, iMody, iModz, &       ! integer,  phys
                         ModAmp,  &                                   ! real,     phys
                         restart, &                                   ! integer,  num
                         pi, iunit                                    ! constants

  implicit none

!-- interface input/output --

  complex, intent(out) :: dfu(jv, jz, jy, jx, jDir)
  complex, intent(out) :: uxMat(jz, jy, jx, jDir)
  complex, intent(out) :: uyMat(jz, jy, jx, jDir)
  real,    intent(out) :: f0vec(jv)
  complex, intent(out) :: finit(jv,jz,jDir)
  integer, intent(out) :: it


!-- local variables ---

  integer              :: ix, iy, iz, iv, iDir
  real                 :: MoDinZ(jz), MoDinX(jx), MoDinY(jy) 
  complex              :: dummy3D(jz,jy,jx)


!-- create vector "f0vec" --


  f0vec(:) = ( 1/sqrt(2*pi) )*exp(-0.5*(vvec(:) + vPhase*iWaveFrame)**2)
                 


!-- read restart or create new --

  if (restart) then

       write(msg,*) char(10), ' Reading restart', char(10)
       call mpi_wrap_msg

       call read_bin(dfu, uxMat, uyMat, it)        


  else
 
       write(msg,*) char(10), ' Creating initial data', char(10)
       call mpi_wrap_msg

       it = 0
       dfu=0.0

       do iDir = 1, jDir
          uxMat(:,:,:,iDir) = ux0(iDir)
          uyMat(:,:,:,iDir) = uy0(iDir)
       enddo
 
       if ( iModx ) then
          MoDinX(:) = cos( dkx * xvec(:) )
       else
          MoDinX(:) = 0
       end if

       if ( iMody ) then
          MoDinY(:) = cos( dky * yvec(:) )
       else
          MoDinY(:) = 0
       end if

       if ( iModz ) then         
          MoDinZ(:) = cos(  kz * zvec(:) )
       else
          MoDinZ(:) = 0
       end if


       do ix=1,jx                                                      
          do iy=1,jy
             do iz=1,jz
                do iv=1,jv
                   dfu(iv,iz,iy,ix,:) =  &
                        ModAmp * f0vec(iv) * ( MoDinZ(iz) + MoDinX(ix) + MoDinY(iy))
                enddo
             enddo
          enddo
       enddo

  end if

!-- save initial (1D, equilibrium) dfu.

  do iDir = 1,jDir
     do iv=1,jv

        dummy3D(:,:,:) = dfu(iv,:,:,:,iDir)

        call fft_forward_x(dummy3D)
        call fft_forward_y(dummy3D)

        fInit(iv,:,iDir) = dummy3D(:,1,1)

     end do
  end do



  !fInit(:,:,:)=dfu(:, :, 1, 1, :)      

end subroutine set_ic


!-------------------------------------------------------------------------------------
subroutine set_ic_lin( dfu, uxMat, Fbgk, it)
!-------------------------------------------------------------------------------------

  use dimensions
  use mpi_wrap
  use fft_3d
 
  use parameters, only : ux0, pi, restart, runname, kz, ModAmp, iSingleMode, iunit

  implicit none


  complex, intent(out) :: dfu(jv, jz, jy, jx, jDir)
  complex, intent(out) :: uxMat(jz, jy, jx, jDir)
  real,    intent(in)  :: Fbgk(jv,jz)
  integer, intent(out) :: it

  integer              :: ix, iz, iv, iDir, filenum
  character*7 :: file_ext

  real        :: epsilon
  complex     :: noise(jz,jy,jx)
  real        :: noise1(jz,jy,jx)
  real        :: noise2(jz,jy,jx)
  complex     :: uyMat(0,0,0,0)

!-- read restart or create new --

  if (restart) then

       write(msg,*) char(10), ' Reading restart', char(10)
       call mpi_wrap_msg

       call read_bin(dfu, uxMat, uyMat, it)        

       return

  end if
 
       write(msg,*) char(10), ' Creating initial data', char(10)
       call mpi_wrap_msg

       it = 0
       dfu=0.0
       uxMat=0.0 

       !call random_number(dfu)
       !dfu(:,:,:,:,:) = (dfu(:,:,:,:,:)*2 - 1)*ModAmp
       !
       !do iDir=1,jDir
       !   do iv=1,jv
       !      noise=dfu(iv,:,:,:,idir);
       !      call fft_forward_z(noise);
       !      noise(1,:,:)=0;
       !      call fft_backward_z(noise);
       !      dfu(iv,:,:,:,idir)=noise;
       !   enddo
       !enddo

       !-- multi-mode
       !
       !call random_number(uxMat)
       !do iDir = 1, jDir
       !   uxMat(:,:,:,iDir) = (uxMat(:,:,:,iDir)*2 - 1)*ModAmp + ux0(iDir)
       !enddo


       !-- single-mode

       noise = 0;
       if (iSingleMode) then
          noise(2,1,1) = 1;
       else
          noise(2,1,2) = 1;
       end if

       call fft_backward_x(noise);
       call fft_backward_y(noise);
       call fft_backward_z(noise);

       do iDir = 1, jDir
          uxMat(:,:,:,iDir) = noise*ModAmp
       enddo

       
       !do iDir = 1, jDir
       !   call random_number(noise1)
       !   call random_number(noise2)
       !   noise = (noise1-0.5) + iunit * (noise2-0.5)
       !   call fft_backward_z(noise);       
       !   uxMat(:,:,:,iDir) = noise(:,:,:) * ux0(iDir) * ModAmp
       !enddo


       !-- single kx,kz mode

       !do iz = 1,jz
       !   do iv = 1,jv
       !      dfu(iv,iz,:,1,:) =   0.01 * cos(kz*zvec(iz)) * Fbgk(iv,iz)
       !      dfu(iv,iz,:,3,:) =  -0.01 * cos(kz*zvec(iz)) * Fbgk(iv,iz)
       !   end do
       !end do


end subroutine set_ic_lin


!-------------------------------------------------------------------------------------
subroutine set_ic_nlvp(dfu, Fbg, it)    ! outputs
!-------------------------------------------------------------------------------------

  use dimensions
  use mpi_wrap
  use fft_wrap 

  use parameters, only : restart, runname, ModAmp, iLinVlasov, kz

  implicit none


  complex, intent(out) :: dfu(jv, jz, jy, jx, jDir)
  real,    intent(in ) :: Fbg(jv)
  integer, intent(out) :: it

  integer              :: ix, iy, iz, iv, iDir, filenum
  integer              :: jvc

  character*7          :: file_ext
  complex              :: fv(jv), fz(jz)
  complex              :: uxMat(0,0,0,0)
  complex              :: uyMat(0,0,0,0)

!-- read restart or create new --

  if (restart) then

       write(msg,*) char(10), ' Reading restart', char(10)
       call mpi_wrap_msg

       call read_bin(dfu, uxMat, uyMat, it)        

       return

  endif


 
  write(msg,*) char(10), ' Creating initial data', char(10)
  call mpi_wrap_msg

  it = 0
  dfu=0.0

  !call random_number(dfu)
  !dfu(:,:,:,:,:) = dfu(:,:,:,:,:)*2 - 1
  
  iDir=1
  ix=1
  iy=1
  
  !-- filteing in z
  !
  !do iv=1,jv
  !   fz=dfu(iv,:,iy,ix,idir)
  !   call dfftw_execute_dft(plan_zforward, fz, fz)
  !   fz = fz/float(jz)
  !   fz(1)=0
  !   fz(3:jz-1)=0                 ! single kz-mode
  !   call dfftw_execute_dft(plan_zbackward, fz, fz)
  !   dfu(iv,:,iy,ix,idir)=fz;
  !enddo


  !-- filteing in v
  !
  !jvc=jv/4        
  !do iz=1,jz
  !   fv=dfu(:,iz,iy,ix,idir)
  !   call dfftw_execute_dft(plan_vforward, fv, fv)
  !   fv = fv/float(jv)
  !   fv(2+jvc : jv-jvc)=0         ! keep only "jv < jvc" modes
  !   call dfftw_execute_dft(plan_vbackward, fv, fv)
  !   dfu(:,iz,iy,ix,idir)=fv;
  !enddo


  !-- normalizing to ModAmp
  !
  !dfu = dfu * ModAmp / maxval(abs(dfu))
 

  !-- single kz mode

  do iz = 1,jz
     dfu(:,iz,:,:,:) = ModAmp * cos(kz*zvec(iz))
  end do

  !-- modulated in z

  if (iLinVlasov) then

     do iv = 1,jv
        dfu(iv,:,:,:,:) = dfu(iv,:,:,:,:) * Fbg(iv)
     end do

  else

     do iv = 1,jv
        dfu(iv,:,:,:,:) = (dfu(iv,:,:,:,:) + 1.0) * Fbg(iv)
     end do

  end if

end subroutine set_ic_nlvp


!-------------------------------------------------------------------------------------

